Astronomy & Astrophysics manuscript no. frlfinal 


© ESO 2008 


June 10, 2008 





Formation of dynamical structures in relativistic 

jets: tlie FRI case 

p. Rossi \ A. Mignone''^, G. Bodo\ S. Massaglia^, and A. Ferrari^'^ 

'66 

o 

■ ' INAF/Osservatorio Astronomico di Torino, Strada Osservatorio 20, 10025 Pino Torinese, Italy 



(N 



X 



^ Dipartimento di Fisica Generale, Universita degli Studi di Torino, via P. Giuria 1, 10125 Torino, 
Italy 

' Department of Astronomy and Astrophysics, University of Chicago, 5640 S. Ellis, Chicago, IL 
60637 



Received ; accepted 

ABSTRACT 

o 

■ Context. Strong observational evidence indicates that all extragalactic jets associated with AGNs 

, move at relativistic speed up to 100 pc - 1 kpc scales from the nucleus. At larger distances, re- 

flecting the Fanaroff-Riley radio source classification, we observe an abrupt deceleration in FR-I 
^ ^ jets while relativistic motions persist up to Mpc scale in FR-II. Moreover, VLBI observations of 

QQ ' some object like B2 1 144+35, Mrk501 and M87 show limb brightening of the jet radio emission 



at the parsec scale. This effect is interpreted kinematically as due to the presence of a deboosted 
central spine at high Lorentz factor and of a weakly relativistic external layer. 
• Aims. In this paper we investigate whether these effects can be interpreted by a breaking of the 

. collimated flow by external medium entrainment favored by shear instabilities, namely Kelvin- 

, Helmholtz instabilities. We examine in details the physical conditions under which significant 

deceleration of a relativistic flow is produced. 

Methods. We investigate the phenomenon by means of high-resolution three-dimensional rela- 



$H ' tivistic hydrodynamic simulations using the PLUTO code for computational astrophysics. 

Results. We find that the parameter of utmost importance in determining the instability evolution 
and the entrainment properties is the ambient/jet density contrast. We show that lighter jets suffer 
stronger slowing down in the external layer than in the central part and conserve a central spine 
at high Lorentz factor. 

Conclusions. Our model is verified by constructing synthetic emission maps from the numeri- 
cal simulations that compare reasonably well with VLBI observations of the inner part of FR-I 
sources. 

Key words. Radio jets - Numerical simulations - Relativistic flows - Shear instabilities 

1. Introduction 

Extragalactic radio s ources are traditionally divided in two morphological classes according to 



their intrinsic power jFanaroff &Rilevll974 : low luminosity sources (Fanaroff-Riley type I, FR-I) 
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are brighter close to the nucleus of the parent galaxy and their jets become dimmer with distance, 
while high power sources (FanarofF-Riley type II, FR-II) show the maximum brightness in the hot 
spots at the jet termination. The different morphology is generally accepted to reflect a difference in 
how the jet energy is dissipated during propagation in the extragalactic medium and produces the 
observed radiation. For FR-I sources, it was quickly accepted that entrainrnent and deceleration of 



the jet must play an imp ortant role in shaping their morphology jBickneU 



1984 



1986; 



1996; 



Komissarov 



De Young 



19941) . while in FR-II sources energy and momentum are transported without 



losses to the front working surface. 

More recently a large body of evidence has accumulated showing that jets are relativistic at 
their base, not only in FR-II radio sources b ut also in FR-I. Superlu minal motions are observed 



on milliarcsecon d scales in several FR-I jets ([Giovannini et al 



M87 and Cen A jBiretta et al 



1995 



Hardcastle et al. 



2001) and on arcsecond scales in 



20031) and finally, on small scales, there are 



also observations of one-sidedness and brightriess as ymmetry between jet and counter-jet, due most 
likely to Doppler boosting effects jLaing et al .I1999I) . FR-I sources are also thought to be the parent 
population of BL Lac objects, for which the presence of relativistic velocities on parsec scales is 



well established (lUrrv & Padovani 



19951) . The Lorentz factors of the jet bulk motion at sub-pc 



scales cannot be deduced direc tly from the observations, but are inferred from assumption on the 
physical emission mechanism. [Harris & Krawczvnskil ( l2006h in their review i ndicate that FR-II 
radiogalaxies have jets with bulk Lorentz factors between 5 and 40 (see also 



Urrv & Padovani 



I995I) and that FR-I jets have less constrained values somewhat lower than those of the FR-II jets. 

( 1200 lb using data on proper motions and brightness ratio between jet and counter 



Giovannini et al 



jet in a sample of radiogalaxies do not find any systematic diffe rence between low and high power 

koom derive 



Celotti & GhiseUini 



radio sources, and the values they derive are between 3 and 1Q.[ 
jet physical properties modelling the spectral energy distribution in a sample of blazars and find 
no difference between BLLac's (associated to low power radio sources in the unified model) and 
radio-quasars (associated to high power radio sources), with an average value of about 15. 

On the other hand, relativistic motions from the inner regions all the way to larger scales, with 
Lorentz factors of about ten, appear to be present i n powerful FR-II jets, as indicated by Chandra 



discovery of bright X-ray emission at kpc scales dTavecchio et al. 



20061) : instead sub-relativistic velocities are found at kpc scales JBicknell 



2004 



Harris & Krawczvnski 



19941) in low power 



radio-sourc es, and the decrease of brightness asymmetry along the jet s uggests that a deceleration 



must occur dLaing et al.ll 19991; iLaing & Bridlj2002 ; 



Canvin et alj|2005h . 



Both these morphological and kinematical data are indicative of the interaction between col 
limated outflows and the surrounding medium. In particular jet decelerati on can be obtained b y 
redistributing the bulk momentum through some form of mass entrainment ( Bicknelll994 



1995). 



Velocity she ar instabilities are the most likely triggering mechanisms of entrainment dDe Young 



1996 . 



ZooJ.Q In fact the nonlinear development of velocity shear or Kelvin-Helmholtz instabili- 



ties leads to an exc hange of mass, momentum and energy at contact discontinuities between fluids 



in relative motion dBodo et al 



1994 



1995 



19981) : understanding the details of this process is es- 



sential for modeling the jet dynamics and for giving clues on the determination of the jet physical 



An alternat i ve and mo st likely complementa ry view for the origin of the entrainment has been presented 



Komissarov 



19941) and 



Bowman et al. 



by 



mass lost by stars within the jet volume. 



1996I) . who consider the possibility of entrainment by injection of 
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parameters. One of the key features of this form of interaction is the formation of a mixing layer at 
the interface between jet and surrounding medium, where the external material is entrained and ac- 
celerated at the expenses of the jet momentum. This process leads to the formation of a transverse 
velocity profile where the internal layers feel less the effects of the interaction with the ambient 
medium and keep an higher velocity, while the external layers are more decelerated. The presence 
of such transverse velocity structure has been already suggested for explaining some of the obser- 



1990; 



vational pro perties of radio sources su ch as their magnetic field configuration jKomissarov 
Laing|[l993 ), limb brightening effects ( Giroletti et al.ll2004 ) and to overcome problems in unifying 



radiogalaxies with BL Lac objects dChiaberge et al 

For studying in detail the instability evolution and the subsequent entrainment process one has 
to resort to three-dimensional numerical simulations, since the mechanisms at the base of t he en- 
trai nment are inherent ly three-dimensional, as shown in preliminary analyses bv lBodo et al.l (120031) 
and kossi et all dlOOJ). 

In this paper we present results of high-resolution hydrodynamic simulations in which we fol- 
low the evolution of a perturbed relativistic jet as it propagates in a homogeneus stationary ambient 
medium. We assume that the jet is in pressure equilibrium with the outside medium, although this 
choice is not critical for the final results. Perturbations grow as a consequence of the velocity shear 
instabilities and lead to entrainment of external medium and to jet deceleration. The main questions 
we address are the dependence of the entrainment and deceleration processes on the jet physical 
parameters and the kind of structure that the jet acquires as the result of these processes. 

The plan of the paper is the following: in ^we present the numerical setup adopted for cal- 
culations and the parameter space covered by the simulations; in ^we discuss the results of our 
simulations focusing our attention on the dependence of the efficiency of deceleration on the physi- 
cal parameters; in ^we analyze in more detail the entrainment properties for the case that appears 
more successful in decelerating the jet. Preliminary comparisons of our simulations with FR-I 
source of different morphologies are reported in ^and conclusions of our study are drawn in ^ 



2. Numerical setup 

Numerical simulations are carried out by solving the equations of particle number and energy- 
momentum conservation. Referring to the observer's reference frame, where the fluid moves with 
velocity Vk (in units of the speed of light c) with respect to the coordinate axis k - x,y,z and 
assuming a flat metric, the conservation laws take the differential form: 





py 




pyvi 






d 
di 


wy^vk 

wy^ - p 

. pyf , 


+ y ^ 

^ 4^ dxi 


wy^VkVi + p6ki 

pyfvi 


= 0, 


(1) 



where p, w, p and y denote, respectively, the rest mass density, enthalpy, gas pressure and Lorentz 
factor. The jet and external material are distinguished using a passive tracer, /, set equal to unity 
for the injected jet material and equal to zero for the ambient medium. The system of equations 
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( B is completed by specifying an equation of state relating w, p and p. Following iMignone et al 
(120051) . we adopt the following prescription: 



:^p+ yj^p^ (2) 

which closely reproduces the thermodynamics of the Synge gas for a single-specie relativistic 
perfe ct fluid offering, at the same ti me, considerable ease of implementation and reduced numerical 
cost ( Mignone & McKinnevlboO? ). 

Simulations are carried out on a Cartesian domain with coordinates in the range x G 
[-Lji/2, Lx/2], y € [0, L^] and z e [-L^/2, L,/2] (lengths are expressed in units of the jet radius, y is 
the direction of jet propagation). At f = 0, the domain is filled with a perfect fluid of uniform density 
and pressure values representing the external ambient medium initially at rest. A constant cylindri- 



cal inflow is initialized in the small region y < I, Vx^ + z} < 1 with velocity along the y coordinate 
and is constantly fed into the domain through the lower y boundary. Jet inflows are represented by 
a three-parameter family distinguished by the beam Lorentz factor yb, Mach number M - Vh/cs 
and the ambient-to-beam density contrast rj as seen in the observer's frame. Thus, by normalizing 
the beam density to phyh = 1, the ambient medium is given by pa = rj/jb- The beam is initially 
pressure-matched and pb is recovered from the definition of the Mach number. This assumption is 
the most used in jet simulations and in this first study we maintain it, however it can be questiond, 
since it is not clear what could be the mechanism that brings pressure equilibrium between the je t 



Komissarov & Falld d 19981) 



and the surrounding medium. An alternative approach is suggested by[ 
and could be worth investigating in future studies. In our simulations we fixed the jet Lorentz fac- 
tor equal to 10 at the jet inlet, a representative value in the range derived from the observations 
discussed in the Introduction. For the Mach number we choose the values of 3 and 30; however for 



relativistic fluids it is actually m ore appropriate to consider the relativistic Mach number jKonigl 



1980; 



Komissarov & Falle 



19981) defined as Mr - ybVbljsCs. where js - 1 / -^/T^^cfTQ ■ It is in 



fact Mr that has the same physical interpretation as the usual non-relativistic Mach number since it 
reflects the propagation of sound waves and the value of Mach angle. Finally for the density ratio 
we choose the values 10^ and lO''. The complete set of parameters for our simulations is given in 
Table 1. 

At the jet inlet, transverse velocities are perturbed introducing pinching, helical and fluting 
modes with corresponding azimuthal (0 = tan"'(z/x)) wave numbers m = 0, 1, 2: 

(Vx, = ^ X X ^'""^ '^'^ ^'^ ^'^"^ ' 

m=0 /=1 

where high (/ - 1, ... ,4) and low (/ - 5, . . . , 8) frequency modes are given, respectively, by 
w/ - Cj(l/2, 1,2,3) and w; = Cj(0. 03, 0.06, 0.12, 0.25). The phase shifts bi are randomly chosen. 
The amplitude A of the perturbation corresponds to a fractional change in the bulk Lorentz factor, 
Jb ^ 7^(1 + e), yielding 

A^-^^^^^y (4) 

Jhd + e) 

We typically set e - 0.05. Outside the inlet region we impose symmetric boundary conditions 
(emulating the presence of a counter-jet), whereas the flow can freely leave the domain throughout 
the remaining boundaries. 
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Yb 


M 


M,. 


'1 


pts/bcam 




A'r X Ny X N. 


A 


10 


3 


28.3 


102 


20 


50 X 150x50 


324 X 1200 X 324 


B 


10 


3 


28.3 


lO'* 


20 


60 X 75 X 60 


344 X 600 X 344 


C 


10 


3 


28.3 


lO'* 


12 


50 X 75 X 50 


172 X 300 X 172 


D 


10 


30 


300 


10"* 


20 


50 X 150 X 50 


324 X 1050 X 324 


E 


10 


30 


300 


10^ 


12 


24 X 200 X 24 


144 X 560 X 144 



Table 1. Parameter set used in the numerical simulation model, the second column refers to Lorentz 
factor jti, the third to the Mach number, the fourth to the relativistic Mach number, the fith to the 
ratio of proper densities, the sixth to number of mesh points on jet radius, the seventh the physical 
domain in units of jet radii and the last one the numerical domain 



Explicit numerical integration of the equat ions dB is achieved u sing the relativistic hydrody- 
namics module available in the PLUTO code (IMignone et al .112007b . PLUTO is a Godunov-type, 
multi-physics code providing a variety of computational algorithms for the numerical integration 
of hyperbolic conservation laws in multiple spatial dimensions. For the present application, we em- 
pl oy the dimensionally s pht version of the relativistic Piecewise Parabolic Method (PPM) presented 



Mignone et al.l ( l2005b . which has an overall 2"'* order accuracy in both space and time. 



The physical domain is covered by NxXNyXN^ computational zones, not necessarily uniformly 
spaced. For domains with a large physical size, we employ a uniform grid resolution around the 
beam (typically for |z| < 2) and a geometrically stretched grid elsewhere. 



3. Results 



The overal l char acteristics of relativistic jet propagation have been originally discussed by 
(Il997h who covered a wide portion of the parameter space wi th 2D simulations. Three- 



Marti et al 



dimensional resu lts for one of the cases studied by 



Marti et al. 



Alov et al 



(119971) have been presented by 



(119991) . Their simulations generalize to the relativistic case typical results of the non- 
relativistic regime, namely: (i) hghter jets exhibit fatter, nearly spherical cocoons (see e.g. Krause & 
Camenzind 2003) as a consequence of the reduced head propagation velocity and bulk momentum 
of the beam; (ii) heavier jets advance at higher velocities generating cylindrical, thinner cocoons; 
(iii) high Mach numbers have the same effect as low density contrasts in prod ucing elongated 



(spea r head) cocoons, wh ile low Mach numbers contribute to form fatter cocoons (IMassaglia et al 



199a) . 



Marti et al 



(Il997h show that a similar behavior is found for relativistic jets, the main differ- 
ence being that, increasing the Lorentz factor, jets, even at low density, tend to form cocoons that 
resemble those of the heavy non-relativistic case due to their increased inertia. In addition in rela- 
tivistic fluid dynamics, different Mach numbers correspond to different jet internal energies rather 
than beam velocities, thus being the density contrast t] the main parameter defining the shape of the 
cocoon. 

Our results, as we shall see, confirm this global trend; however, in this scenario, we intend 
to focus our attention on the entrainment properties as a possible source of the jet deceleration. 
The entrainment process takes place from the interaction between the jet beam and the cocoon, 
promoted by the development of Kelvin-Helmholtz instabilities (KHl henceforth) at the beam in- 
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Fig. 1. Volume rendering of the tracer distribution for case A at / = 240 (top panel), B at / = 760 (central 
panel) and E at ? = 150 (bottom panel). Dark colors are due to the opacity of the material in the rendering 
procedure 



terface. The entrained material is composed by jet backflowing material mixed with the shocked 
ambient medium through the contact discontinuity. The process is therefore quite complex and 
determined by many different factors, the behavior of KHI being of course one of the most rele- 



vant. The linear analysis of the KHI for relativistic flows (lHardeell987 



ng of course one of t 
:ell987[lHardee et al.ll998t 



Hardee 



20001) shows that the growth rate of unstable modes depends on the Mach number, density ratio 



between jet and external medium and Lorentz factor. In this respect, we intend to investigate for 
which set of parameters entrainment occurs more efficiently. 

Obviously the numerical resolution is a very important factor in simulating the mixing pro- 
cesses. Therefore we perform our study with very high-resolution 3D computations in the hydro- 
dynamic regime. The presence of magnetic field may change the global structure, but as long as it 
is not too strong it is not believed to change the diffusion processes. Given the high computational 
cost of 3D simulations, we had to limit our study to a fixed value of the Lorentz factor varying 
only the Mach number and density ratio. Table 1 lists the complete set of parameters adopted in 
our numerical simulations, together with the mesh resolution. With the choice of these elements we 
can cover a complete set of combinations of Mach numbers and density ratios (cases A, B, D and 
E). Moreover, case C represents a lower resolution version of case B, and a comparison between 
the two can provide indications on the role played by numerical resolution. In agreement with ob- 
servations, we consider only jets with density much lower than the external medium; on the other 
hand high density jet propagate quasi ballistic ally, do not show instabilities and mixing 
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As we shall discuss below, the main parameter governing the jet behavior and the entrainment 
properties results to be the density ratio rj, while the Mach number plays a minor role even though it 
contributes to the determination of the KHI modes linear growth rates. Our discussion will therefore 
refer to cases A and B, as representatives of low and high i] values. Our results show that case D 
has a behavior very similar to case B, for this reason we will not discuss it in detail. Case E will be 
instead properly discussed in order to show what can be the effects of a higher Mach number. 

Fig. [T] shows the volume rendering of the tracer distributions for case A (top panel), B (mid 
panel) and E (bottom panel) when the evolutions have reached fully nonlinear stages. Jet material 
is dark brown, external medium material white, level of mixing corresponds to the scale of orange. 
From these figures one can gain a first qualitative view about the role of the parameters on the 
evolution and structure of the jet: in cases A and E the jets propagate straight with moderate and 
low dispersion of jet particles, whereas in case B the jet has, at least partially, lost its collimation 
and the particle dispersion is considerably larger. 

This is more clearly represented by the 2-D cuts in the x, y plane (at z = 0) of the density and 
Lorentz factor distributions shown in Fig. |2] From the density (left) panels one can note that in 
both cases A and E the jet seems to be weakly affected by the perturbation growth and entrainment. 
This is not the situation for case B, where the beam structure is heavily modified by the growth of 
disturbances beyond ~ 20 jet radii. The distributions of the Lorentz factor (right panels) indicate 
again that the perturbation slightly affects the system in case A after ~ 100 radii, leaves the jet 
almost unchanged in case E, and strongly influences the propagation after ~ 20 radii for case B, 
where the maximum value of y has reduced approximately to half of its initial value. However, 
even in case B a well coUimated high velocity component along the axis still displays. 

A quantitative estimate of the jet deceleration can be obtained by plotting the Lorentz factor 
as a function of the longitudinal coordinate y. Figure [3] shows the maximum value of y at constant 
y-planes together with its volume average, defined as 

/ ygij) dxdz 

Tav = -7 , (5) 

J giy) dxdz 

where giy) is a filter function to select the relativistic flows: 
f 1 for y > 2 , 

g{y) = (6) 
\Q for y < 2 . 

For case B one can see that the deceleration occurs both in ymax (the flow velocity, although still 
relativistic, shows a strong decrease from its initial value) and in y-^y, which indicates a global 
effect, whereas in cases A and E the central part of the jet continues to be almost unperturbed and 
only thin external layers are decelerated. 

As previously mentioned, the interaction between jet and surrounding medium involves differ- 
ent phenomena, reflecting on the characteristics of the deceleration process: the growth rate of the 
perturbations, the type of perturbations that dominates the jet structure, the possibility of mixing 
through the backflowing material and the density of the ambient medium relative to the jet material. 
Long wavelength modes will tend to produce global jet deformations like jet wiggling more than 
mixing, while modes at shorter scale may be more efficient in promoting the mixing process. On 
the other hand, it is clear that the mixing with a denser environment is more effective for the jet 
deceleration. One of the differences between the two cases that we are considering is, in fact, the 
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Fig. 2. Two-dimensional longitudinal cut in the xy plane of the density distribution for case A at f = 240 
(top-left panel), case B at f = 760 (central-left panel) and case E at f = 150 (bottom-left panel), and Lorentz y 
distribution in the central part of the domain, not in scale (the corresponding right panels). 



type of modes that dominate the jet evolution: while case E exhibits longer wavelength and lower 
m modes with a highly prominent helical distortion, in case B we observe shorte r wavelength and 
higher m modes. This behavior agrees with the linear analysis of lHarded ( Il987h that predicts that 
wavelengths of the resonant mode increases with the Mach number as /I oc jM. An impression of 
the characteristic scales of the dominant modes can be obtained from Fig. H) where 3D contour 
images of the Lorentz factor are shown. The central panel, case B, shows that small scale structures 
dominate the structure all along the jet, while case E (right panel) shows the dominance of the 
m - 2 helical mode. The left panel (case A), denotes a less clear cut situation possibly due to the 
effect of the superposition of different modes, with the presence of structures of scale smaller than 
in case E. 

So far we have compared our three simulations at a given time. In Figs. |5] |6l |7] [8l |9] and 
[TOl we consider as well the temporal evolution of pressure, density and velocity, comparing the 
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Fig. 3. Plots of tiie maximum value (solid line) and of an average value (dashed line) of y as functions of the 
longitudinal coordinate y along the jet. The left panel refer to case A, the center panel to case B and the right 
panel to case E. 






Fig. 4. Image of the Lorentz factor distribution of the contour for case A (left panel, y = 5), case B (central 
panel, 7 = 3) and case E (right panel, y = 5). 



jet structures at three different times when they reach approximately equal lengths. In Fig. |5] 
we show cuts in the plane of the pressure distribution with superimposed contours of the 
Loren tz factor. The figure displa ys the formation of the typical series of regular biconical struc- 
tures (IKomissarov & Fallelll998h produced by the overpressured cocoons squeezing the jet. Case 
B shows the highest pressure and has shocks of smaller obliquity. From the contour levels of the 
Lorentz factor we see once more that cases A and E keeps the initial value of 10 up to the jet head. 
Fig. |6] shows the jets when they have reached a length of 120 radii; here we display only cases A 
andE since, for case B, the jet head velocity is very low and the simulation was ended at a jet length 
of 60 radii. In the figure, we start seeing differences between case A, more perturbed, and case E, 
that shows almost no sign of interaction between the jet beam and the cocoon. Differences between 
cases A and E were aheady evident in Fig. |4] that showed the presence of small scale perturbations 
more favourable to mixing in case A and of only a long wavelength hehcal distortion in case E. 

On the contrary, for case B that shows the strongest deceleration, a first decrease of the Lorentz 
factor occurs at the first biconical shock (see also Figs. |2] and |3]l, the reacceleration following the 
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Fig. 5. Transverse cuts in the x - y plane of the pressure distribution with superimposed the Lorentz factor 
contours. The left panels refer to case A, the middle panels to case B and the right panels to case E. The upper 
panels show the jet when it has reached a length of 30rj while, in the lower panels, the jet length is 60 rj. 

shock however does not lead to a complete recovery of the initial value because the jet beam starts 
to interact with the cocoon and to lose momentum. We can note that, in the relativistic case, an 
increase of the jet pressure to resist squeezing by the overpressured cocoon would result in an 
increase of the jet energy flux (due to the relativistic contribution of pressure to inertia) and in a 
further increase of the cocoon pressure. Thus, in the relativistic regime, the formation of biconical 
shocks appears to be even more unavoidable than in the nonrelatistic regime. 

The properties of the entrainment process can be better understood by looking at Figs ITlfTOl 
In Fig. |7] and [8] we show cuts in the ix,y) plane of the density distribution of the external medium 
(derived from the tracer distribution) for the three cases considered, at different times. This figure 
shows the ability of the external medium to penetrate the contact discontinuity and mix with the jet 
beam. We see that in B (central panel) the jet is surrounded by denser (green-red) external material 
in most of its length, while in cases A (left panels) and E (right panels) the medium around the 
jet is always lighter (blue-green). Comparing the different evolutionary stages (see Figs [T] and [8]l, 
we also observe a slight decrease in density. The effect of the jet on the speed of entrained matter 
is shown in Figs. |9] and (TO] where we display cuts of the distribution of the external medium 
velocity component along the jet axis. We note that the volume of the accelerated medium (red) 
is considerably higher in case B than in cases A and E, consistently with a substantially stronger 
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Fig. 6. Transverse cuts in the x-y plane of the pressure distribution of the external material with superimposed 
the Lorentz factor contours. The left panel refers to case A, the and the right panel to case E. The jet length is 
120r, 



U = V7 - 10- \f = 3 Ti - 10^ ^f = 30 77 - in" 




.■if; ■<; o r-i 20 -ao ;»i 20 <(i o vd :iv, ■-:>.) :•:(! 20 ui w jo 2r, no 



T J' 

Fig. 7. Transverse cuts in the x-y plane of the density distribution of the external material. The left panels 
refer to case A, the middle panels to case B and the right panels to case E. The upper panels show the jet when 
it has reached a length of 30rj while, in the lower panels, the jet length is 60 rj. 
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Fig. 8. Transverse cuts in the x — y plane of the density distribution of the external material. The left panel 
refers to case A, the and the right panel to case E. The jet length is 120?; 

instability. Entrainment in case B starts to be effective already at a length of 30 radii, when it is 
almost absent in the other cases. At a length of 60 radii, apart from the dominant case B, we start 
observing differences between case A and case E with a stronger entrainment effect of the first. 
This difference becomes even more prominent at a length of 120 jet radii (see Fig.fTOl). From this 
figures we can also observe the formation of the backflow (blue region), more prominent in case B 
and almost absent in cases A and E. 

4. Entrainment properties 

We now discuss in detail the entrainment properties of the jet of case B (M - 3,t] - 10"^). We begin 
examining the distribution in Lorentz factor and velocity of the jet and external material at f = 600 
(at the end of the simulation). 

In Fig.[TT]we show the distribution of the jet mass fraction that is moving at a certain value of 
yjS at a given time. We use the four-velocity (instead of the Lorentz factor) to avoid compression of 
the scale close to y = 1, i.e. at low velocities. The color scale on the right gives the corresponding 
value of y/S for each color. As an illustrative example, one can see that at y = 20, ~ 20% of the 
jet mass (per unit length) moves at 4 < y/? < 5, ~ 40% moves at 3 < y/3 < 4, ~ 17% moves at 
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Fig. 9. Transverse cuts in the x - y plane of the v,, distribution of the external medium. The left panels refer 
to case A, the middle panels to case B and the right panels to case E. The upper panels show the jet when it 
has reached a length of 30ry while, in the lower panels, the jet length is 60 rj. 

2 < j/3 < 3, ~ 5% moves at 1 < yj3 < 2 and finally ~ 18% moves at 0.1 < y/3 < 1. In this image 
one can clearly recognize the effect of reacceleration after internal shocks driven in the beam by 
the high pressure cocoon, particularly close to the inlet region. The central region of the jet never 
decelerate below 7 = 4, while in the external layers of the jet we observe the formation of an 
expanding region of low velocity material {y/3 < 1). 

Additional details are illustrated in Fig. [12] where the jet mass distribution is plotted as a func- 
tion of -y/S at three different values of y for two different instants. The three panels, from left to 
right, refer respectively to y = 12.5, y = 25 and y - 37.5, at f = 400 (dashed line) and t = 600 
(solid line). At the first position {y = 12.5) the differences between the two times is relatively small, 
while at the two other locations we have a significant increase of the material moving at high y. 
The jet in its first part has reached a quasi-steady configuration, while at larger distances from the 
source, after an initial deceleration phase, a well collimated high velocity core appears. Looking at 
the velocity structure we see that, at the beginning (y - 12.5, see the left panel of Fig.[T2li. the jet is 
predominantly at high values of yyS, even though the deceleration mechanism starts to be effective, 
as shown by the small peak at y/3 ~ 0.2. The effects of entrainment become much more evident at 
larger distances (see the two rightmost panels in FigfT2b. where we observe the formation of two 
sharp peaks, one at high y/3 and the other at velocities yP ~ 0.2, with considerably less material at 
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Fig. 10. Transverse cuts in the x-y plane of the v,, distribution of the external medium. The left panel refers 
to case A, the and the right panel to case E. The jet length is llOrj 
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Fig. 11. Distribution of the jet mass fraction moving at a given y/?, when the jet length is 60rj 

intermediate velocities. This implies that the jet structure has two well defined velocity components 
with a steep shear layer between them. 



0^ 4 68/0 Z 4 6 B 10 

7fi yP 



Fig. 12. Jet mass distribution as a function of yj} at the locations y = 12.5, 25, 37.5 at t = 400 (dashed line) 
and t = 600 (solid line) 
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Fig. 13. The same as in Fig.[TT]for tiie external medium and normalized with Mj(0) 

The distribution of external mass, analogous to Fig.[TT|is also plotted , in Fig.[T3] In this case 
all quantities are normalized to the jet mass injected in the domain per unit time. 

The figure show that most of the external material is moving at quite low velocities, i.e. 0.1 < 
/3 < 0.3. A similar analysis on the backflow reveals that most material moves with velocity in the 
range -0.3 < yS < -0.1, with a strongly dominant contribution given by the external medium. 
Higher velocity fluid elements with velocities up to -0.8 can be found in the front regions of the 
jet again with a dominant contribution by the ambient medium. 

We conclude this section with a short evaluation of the effect of numerical resolution. Case C 
has the same parameters as case B, but only 12 points per jet radius, compared to 20 points per 
jet radius of case B. Comparing the results obtained at these two different resolutions, we have 
obtained that the general effect on the entrainment is to increase its efficiency when we increase the 
resolution. A detailed comparison shows that, for example, the average Lorentz factor decreases 
by about 15% in the higher resolution run. This effect can be interpreted observing that the most 
effective modes in term of entrainment are those at a shorter wavelength, that are under-resolved in 
a low resolution run. 



5. Astrophysical implications 



Recently several authors, e.g. jChiaberge et al.ll2000t 



Finer & Edwardsll2004 



Giroletti etal.ll2004) 



in order to explain observational properties of FR-I radio sources and their beamed counterparts 
(BL Lac objects), have proposed that they are produced by jets characterized by a velocity structure 
in which an inner core maintains a highly relativistic velocity and is surrounded by material that has 
been slowed down by the interaction with the ambient medium. A structure of this type has been 
called "spine-layer". The appearance of a jet with a spine-layer configuration is different when 
viewed at different angles. In fact, the two velocity components have different Doppler factors and 
the spine dominates the emission when the jet is observed at small angles with respect to the line 
of sight (BL Lac objects with strong Doppler boosting), while the prevailing contribution at larger 
angles is due to the entrained layer at low Lorentz factors (FR-I radio sources) 

In our calculations, a "spine-layer" velocity structure has been obtained self-consistently as the 
result of a well defined physical process, i.e. the interaction of the outer jet layers with the ambient 
material, driven by jet instabilities. In particular we have found that, in the strongly underdense 
case T] - 10*, the jet acquires a velocity structure in which the inner core maintains a highly 
relativistic velocity and is surrounded by material that has been slowed down by the interaction 
with the ambient medium. Therefore we attempt a comparison of radio maps constructed from 
the simulated jets with observations of FR-I jets. To this purpose, we compute synthetic maps 
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Fig. 14. Synthetic map for the jet of case B with an inclination of 20° with respect to the line of sight. 
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Fig. 15. Synthetic map for the jet of case B with an inclination of 60° with respect to the line of sight. 

by integrating the synchrotron emissivity along the line of sight. For the sake of simplicity we 
assume the emissivity to be proportional to the proper density of the jet material multiplied by the 
appropriate boosting or deboosting factor, i.e. by the quantity 

-\2+a 



eix,y,z) = 



1 



(7) 



r(i -jScose) 

where a is the spectral index of the radio flux (we take a = 0.5). From this calculation we exclude 
the material at rest accumulated at the base of the jet as a result of the reflecting boundary condi- 
tion imposed at y = 0. The presence of this material is an artifact due to the imposed equatorial 
symmetry. 

The synthetic maps for the strongly decelerated jets of case B are plotted in Figs. [T4l and [TS] for 
inclination angles of the jet axis to the line of sight of 20° and 60° respectively. The maps show a 
region with equal projected lengths of 45 jet radii, excluding the jet head, and therefore correspond 
to a longer jet in the small angle case. In the small angle case, the relativistic jet y ~ 4 5 core 
emission dominates. Emission is present along the whole jet and some knots can be recognized, 
the most brilliant being the first one; knots correspond to the oblique shocks in the pressure maps. 
Increasing the inclination angle, the dominant contribution to emission is due to relatively slow 
material, whose fraction in the first part of the jet is very low, so that emission is almost absent. 
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Fig. 16. VLBA radio image at 2 cm of M87 jKovalev et alJl200' 



apart from the first knot that is still visible at lower brightness. The fraction of slow material starts 
to increase at a distance of about 15 radii where the jet appears to have a sudden increase in opening 
angle, due to the formation of a thick layer of slowly moving material which dominates the emission 
at larger distances. We can also observe some limb brightening, explained by the emission coming 
from the slow layer surrounding the deboosted relativistic core. 

We can compare these two synthetic maps with the radio maps of two typical FR-I sources 
shown in Figs. [16] and [T7] The VLBA map of M87 can be compared with our synthetic map at low 
inchnation angle, which agrees with observational estimates. Instead the VLBI map of B2 1 144+35 
can be compared with o ur synthetic map at inclination of 60°, although the estimated inclination 
( Giovannini et al jLoO? ) for this source is 33° ± 7°. One of the possible origins of this discrepancy 
lies in the fact that the Lorentz factor of the slow layer in B2 1 144+35 is estimated to be around 
2.9 (the fast spine has y ~ 15) corresponding to /3 ~ 0.94, while in our simulations (as shown in 
the previous section) the slow material is subrelativistic with /3 ~ 0.3 - 0.4 thus becoming visible 
at larger angles. 



6. Discussion and summary 

In this paper we have presented the 3D nonlinear dynamical evolution of relativistic light jets, as 
a result of a perturbation introduced at the jet inlet. The perturbation grows because of KHI and 
gives rise to a strong interaction of the jet with the external medium with a consequent mixing and 
deceleration. The two main parameters controlling the jet dynamics are the Mach number M and 
the density ratio rj between the ambient medium and the jet. The Lorentz factor has been set equal 
to 10 in all computations. 

We have explored the parameter plane (M, 77), finding that the deceleration becomes more ef- 
ficient increasing 77. A preliminary analysis of the parameter space suggests that only jets with a 
large density ratio (77 > 10^) can undergo appreciable deceleration, while the Mach number does 
not seem to play a fundamental role in this respect. 

We have focused our attention on three extreme cases in the parameter plane, namely: case A 
with M = 3 and 77 = 10^ B with M = 3 and 77 = 10"^ and E with M = 30 and 77 = lOl The com- 




parison of these simulations show in fact that case B (and D) undergoes the strongest deceleration. 
Cases A and E retain a high Lorentz factor spine with propagation velocities essentially unchanged 
from the injection value, while some deceleration has been observed only in the outer layers with 
the formation of a wide transverse velocity structure. In case B we observe the formation of a 
similar pattern although on a much shorter distance and with a significant stronger decrease of the 
maximum Lorentz factor 

The fact that larger values of rj (i.e. lower jet densities) lead to prominent deceleration may 
have direct astrophysical imphcations. Observational data seems to indicate that the jet kinetic 
power a ssociated with FR-I radio sources is, on the average, ~ 10^ times lower than in FR-II radio 
sources JcelottJboOsL On the other hand, there seems to be no difference in the value of the initial 



Lorentz factor in the two classes (iGiovannini et al. 



2001 



Celotti & Ghiseninill2008l) . Since lighter 



jet beams imply reduced jet kinetic powers, our model leaves the density contrast as the most 
likely candidate to account for the discrepancies in the deceleration process efficiency. Using some 
astrophysically relevant units we can rescale our models and come up with a rough estimate for the 
critical value of the jet power f * that separates FRJ from FRJI radiosources: 



p* 
j 



1.3 X 



(8) 
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where we ass umed a jet radius of Ipc , an external density of Icm^^ which is typical at distances 
below IQQpc (IBalmaverde et al.ll2008h . and a critical density ratio separating the two behaviors of 
order 10^. The value o f P* is affected by many un certainities but is in agreement, for ex ample , with 



the estimates given bv lGhisellini & Celottil (1200 lb based on the results by 
gave for P* the value 



Willott et al 



(1 1999b . who 



10" 



'I. 



BH 



ergs 



(9) 



where Mbh is the mass of the central black hole. 

Considering the case with 77 = 10^ we obtained a spine-layer structure similar to that deduced 
from observations. Looking at the synthetic maps produced from the simulations, it is evident that 
many of the salient features are fairly well reproduced. The main difference lies in the terminal 
Lorentz factor of the slow layer, typically smaller than the observational estimates. In view of this 
first promising result we intend to further pursue this investigation, trying to better constrain jet 
parameters and introducing another essential ingredient: the magnetic field. 
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